# ------ Libraries ------ #
library(xergm)
library(ggplot2)
library(superheat)
library(cccd)
set.seed(1912)

setwd("YOUR WORKING DIRECTORY")




# --- Load data --- #

sp_sims <- readRDS("data/word_vectors/rwmds/un.rwmd.l_50d_15x.rds") 
vo_sims <- readRDS("data/comm_detection/candidate_layers/weak/plot_vo_cands.rds")

sp_cands <- readRDS("data/comm_detection/candidate_layers/weak/sp_cands_50d_15x.rds")
vo_cands <- readRDS("data/comm_detection/candidate_layers/weak/vo_cands.rds")

multi_comms_w <- readRDS("data/comm_detection/detected_comms/weak/weak_ties_50d_15x.rds")






# --- Similarity heatmap --- #

t <- 1973

colf <- colorRampPalette(c("gray100", "lightsteelblue4"))

#png("sp_sims.png", width=1000, height=1000)
superheat(sp_sims[[t]], pretty.order.cols = F, pretty.order.rows = F, row.dendrogram = F, col.dendrogram = F,
          left.label = "none", bottom.label = "none", legend = F, heat.pal = colf(20))
dev.off()

#png("vo_sims.png", width=1000, height=1000)
superheat(vo_sims[[t]], pretty.order.cols = F, pretty.order.rows = F, row.dendrogram = F, col.dendrogram = F,
          left.label = "none", bottom.label = "none", legend = F, heat.pal = colf(20))
dev.off()



# --- Network layers --- #

node_colors <- rep("",dim(sp_cands[[t]])[1])
for(i in 1:length(node_colors)){
  cols <- ifelse(rowSums(sp_cands[[t]]) > 0, node_colors[i] <- "lightsteelblue4", node_colors[i] <- "gray100")
}
#png("speech_net_layer.png", width=1000, height=1000)
par(mar=c(.2,.2,.2,.2))
set.seed(123)
plot(network(sp_cands[[t]], directed = F), displayisolates=T, vertex.col=cols, 
     edge.col="gray50", vertex.lwd = 1.8, vertex.border="gray14",
     edge.lwd=5, vertex.cex=2)
dev.off()



node_colors <- rep("",dim(vo_cands[[t]])[1])
for(i in 1:length(node_colors)){
cols <- ifelse(rowSums(vo_cands[[t]]) > 0, node_colors[i] <- "lightsteelblue4", node_colors[i] <- "gray100")
}
#png("vote_net_layer.png", width=1000, height=1000)
par(mar=c(.2,.2,.2,.2))
set.seed(123)
plot(network(vo_cands[[t]], directed = F), displayisolates=T, vertex.col=cols, 
     edge.col="gray50", vertex.lwd = 1.8, vertex.border="gray14",
     edge.lwd=5, vertex.cex=2)
dev.off()


node_colors <- rep("",dim(multi_comms_w[[t]])[1])
for(i in 1:length(node_colors)){
cols <- ifelse(rowSums(multi_comms_w[[t]]) > 0, node_colors[i] <- "#AD6563", node_colors[i] <- "gray100")
}
#png("multi_net.png", width=1000, height=1000)
par(mar=c(.2,.2,.2,.2))
set.seed(123)
plot(network(multi_comms_w[[t]], directed = F), displayisolates=T, vertex.col=cols, 
     edge.col="gray50", vertex.lwd = 1.8, vertex.border="gray14",
     edge.lwd=5, vertex.cex=2)
dev.off()








# --- Weak signal communities over time:

dir <- "data/comm_detection/detected_comms/weak/"
file_list <- list.files(dir, pattern = "weak_membs")


plot_list_w <- list()
for(k in 1:length(file_list)){
  
  multi_comms_w_membs <-  readRDS(paste(dir, file_list[[k]], sep = ""))
  
  plot_df <- data.frame()
  for(i in 1970:1990){
    ncomm <- ncol(multi_comms_w_membs[[i]])
    perc.comm <- length(which(rowSums(multi_comms_w_membs[[i]])>0))/nrow(multi_comms_w_membs[[i]])
    perc.bridges <- length(which(rowSums(multi_comms_w_membs[[i]])>1))/nrow(multi_comms_w_membs[[i]])
    temp <- data.frame(ncomm, perc.comm, perc.bridges)
    plot_df <- rbind(plot_df, temp)
  }
  plot_list_w[[k]] <- plot_df
  
}


plot_df_w <- (plot_list_w[[1]] + plot_list_w[[2]] + plot_list_w[[3]] + plot_list_w[[4]] + plot_list_w[[5]] + plot_list_w[[6]])/length(plot_list_w)
plot_df_w$year <- 1970:1990
plot_df_w$perc.comm <- plot_df_w$perc.comm*100
plot_df_w$perc.bridges <- plot_df_w$perc.bridges*100

ggplot(plot_df_w, aes(x=year, y=ncomm)) +
  geom_smooth(method = "loess", fill = "gray84", colour = "#AD6563", linetype = 2, size = .6) +
  geom_point(shape = 21, stroke = 1, aes(size = perc.comm, fill = perc.bridges)) +
  scale_fill_gradient(low = "#a0b2c9", high = "#404750", 
                      breaks = c(0, round(max(plot_df_w$perc.bridges)/2), round(max(plot_df_w$perc.bridges))),
                      limits = c(0,max(plot_df_w$perc.bridges))) +
  theme_bw() +
  scale_radius(name = waiver(), range = c(2,12)) +
  coord_cartesian(ylim = c(0, 11.8)) +
  guides(fill = guide_colorbar(reverse=T, barwidth = 0.8, ticks = F, 
                               frame.linetype = 1, frame.linewidth = 1.5, frame.colour = "black", 
                               order = 1, label.hjust = .3),
         size = guide_legend(override.aes = list(shape = 1,  alpha = .7))) +
  labs(y="# Communities", x = "Year", fill = "% Bridges", size = "% States in \n a Community") +
  theme(text = element_text(size = 24), legend.title = element_text(size=16), legend.title.align=0.5,
        legend.text = element_text(size=14), 
        axis.title.y = element_text(margin = margin(r=8)),
        axis.title.x = element_text(margin = margin(t=8)),
        panel.grid.major = element_line(colour = "gray84"),
        panel.grid.minor = element_blank())

#ggsave("weak_sig_time.pdf", width = 9, height=5.5)








# --- Strong signal communities over time:

dir <- "data/comm_detection/detected_comms/strong/"
file_list <- list.files(dir, pattern = "strong_membs_masternodes")


plot_list_s <- list()
for(k in 1:length(file_list)){
  
  multi_comms_s_membs <-  readRDS(paste(dir, file_list[[k]], sep = ""))
  
  plot_df <- data.frame()
  for(i in 1970:1990){
    ncomm <- ncol(multi_comms_s_membs[[i]])
    perc.comm <- length(which(rowSums(multi_comms_s_membs[[i]])>0))/nrow(multi_comms_s_membs[[i]])
    perc.bridges <- length(which(rowSums(multi_comms_s_membs[[i]])>1))/nrow(multi_comms_s_membs[[i]])
    temp <- data.frame(ncomm, perc.comm, perc.bridges)
    plot_df <- rbind(plot_df, temp)
  }
  plot_list_s[[k]] <- plot_df
  
}


plot_df_s <- (plot_list_s[[1]] + plot_list_s[[2]] + plot_list_s[[3]])/length(plot_list_s)

plot_df_s$year <- 1970:1990
plot_df_s$perc.comm <- plot_df_s$perc.comm*100
plot_df_s$perc.bridges <- plot_df_s$perc.bridges*100



ggplot(plot_df_s, aes(x=year, y=ncomm)) +
  geom_smooth(method = "loess", fill = "gray84", colour = "#AD6563", linetype = 2, size = .6) +
  geom_point(shape = 21, stroke = 1, aes(size = perc.comm, fill = perc.bridges)) +
  scale_fill_gradient(low = "#a0b2c9", high = "#404750", 
                      breaks = c(0, round(max(plot_df_s$perc.bridges)/2), round(max(plot_df_s$perc.bridges))),
                      limits = c(0,max(plot_df_s$perc.bridges))) +
  theme_bw() +
  scale_radius(name = waiver(), range = c(2,12)) +
  coord_cartesian(ylim = c(0, 8)) +
  guides(fill = guide_colorbar(reverse=T, barwidth = 0.8, ticks = F, 
                               frame.linetype = 1, frame.linewidth = 1.5, frame.colour = "black", 
                               order = 1, label.hjust = .3),
         size = guide_legend(override.aes = list(shape = 1,  alpha = .7))) +
  labs(y="# Communities", x = "Year", fill = "% Bridges", size = "% States in \n a Community") +
  theme(text = element_text(size = 24), legend.title = element_text(size=16), legend.title.align=0.5,
        legend.text = element_text(size=14), 
        axis.title.y = element_text(margin = margin(r=8)),
        axis.title.x = element_text(margin = margin(t=8)),
        panel.grid.major = element_line(colour = "gray84"),
        panel.grid.minor = element_blank())

#ggsave("strong_sig_time.pdf", width = 9, height=5.5)

